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RELIABILITY-BASED FAILURE ANALYSIS OF BRITTLE MATERIALS 


ABSTRACT 


by 

Lynn M. Powers and Louis J. Ghosn 

The reliability of brittle materials under a generalized 
state of stress is analyzed using the Batdorf model. The model 
is modified to include the reduction in shear due to the effect 
of the compressive stress on the microscopic crack faces. The 
combined effect of both surface and volume flaws is included. 
Due to the nature of fracture of brittle materials under compres- 
sive loading, the component is modeled as a series system in 
order to establish bounds on the probability of failure. 

A computer program was written to determine the probability 
of failure employing data from a finite element analysis. The 
analysis showed that for tensile loading a single crack will be 
the cause of total failure but under compressive loading a series 
of microscopic cracks must join together to form a dominant 


crack. 
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CHAPTER I 


INTRODUCTION 


The demand for high -temperature structural materials 
has increased due to the growth in the aerospace, defense 
and energy related industries. Because of the attractive 
physical and mechanical properties of modern ceramics: 
high -temperature strength, lightweight, excellent erosion, 
corrosion and wear resistance, and low thermal conductivity, 
ceramic components are being considered for structural 
applications. Examples of these high -temperature 
applications include: turbine engine components, rocket 
nozzles and nose tips, nuclear fuel pellets and bearings. 

As is the case for other brittle materials, ceramics 
exhibit a large variation in fracture stress which must be 
taken into account in design. This variation in strength 
results from the presence of microscopic random 
imperfections or flaws. Ceramic components contain two 
types of flaws: volume flaws and surface flaws. Volume 
flaws arise from material processing while surface flaws 
arise from grinding and other surface finishing operations. 
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Most of the probabilistic approaches to brittle 
fracture are formulated for tensile failure.^ - However, 
brittle materials react distinctly under different loading 
conditions. These materials tend to be weak in tension, 
but strong in compression. This behavior would suggest that 
two different failure mechanisms govern the fracture of 
brittle materials when in tension or in compression. The 
purpose of this research was to develop a probabilistic 
model to determine the reliability of brittle materials 
accounting for the two mechanisms of fracture and which may 
be applied to any given loading condition, in particular 
contact stress conditions where the stresses are 
predominately compressive. 

A. Reliability Theory 

Reliability is given as the probability of an object 
performing its required function for a specific period of 
time under stated conditions. Vhen considering the failure 
of materials, reliability is measured as the ability of a 
component to sustain load. The statistical nature of 
fracture may be presented in two different theories: the 

O 

weakest link and bundle models. These two concepts are 
illustrated in Fig. 1. 

The weakest link theory makes the analogy between the 
links of a chain and the volume elements of a bulk specimen. 
The strength of the chain is that of its weakest link, the 
strength of the bulk specimen is that of its weakest volume 
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SERIES SYSTEM - WEAKEST LINK MODEL 


PARALLEL SYSTEM - BUNDLE MODEL 



Fig. 1.1 System configuration. 
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element .^* 5 The strength of its weakest volume element is 
characterized by its most severe crack. The weakest link 
theory implies that when one crack fractures the entire 
specimen fails, representative of a series system. 

The alternate concept is that of a bundle model or a 
parallel system. The volume elements within the material 
are composed of links arranged in parallel. Vhen one of the 
links fails, the load is redistributed and the structure may 
survive. Failure is defined when all of the links have 
failed. 

B. Probabilistic Models For Brittle Fracture 

The probabilistic models to be considered are those 
concerned with the instantaneous fracture as a result of 
unstable crack propagation when an initial load is applied. 
Fast fracture of structural components is generally assumed 
to depend on some property of the material from which the 
part was made. It is assumed that brittle materials contain 
microcracks which are uniformly distributed and randomly 
oriented. 

The probabilistic models for brittle fracture are based 
on the weakest link concept, because one crack almost always 
produces total failure in tension. The first of these was 
introduced by Veibull.^’^ Veibull assumed that the 
component of stress normal to the plane of the crack was the 
only one to contribute to its fracture. As a result the 
shape of the crack is irrelevant and crack-crack interaction 
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is not taken into account. The statistical distribution 
function was given as 


a > a. 


u 


* < % 


P r = 1 - exp 


= 0 


L W « 


1.1 


where the probability of failure, P^, is a function of the 
maximum principal tensile stress and three material 
constants: m, <r Q and <r u , the Veibull modulus, scaling factor 
and the threshold stress, respectively. A similar 
distribution for surface flaws may be written by integrating 
over the area instead of the volume. An example of this 
probability distribution is shown in Fig. 2. This theory 
was initially used for materials under uniform tension and 
was extended for non-uniform uniaxial stress states. 

The flaw density distribution characterized by the 
material constants in the Veibull analysis may be difficult 

Q 

to determine. A unique combination of values for the 
three parameters does not exist within "the range used for 
experimental data, 0.05 < P^ < 0.95. For a single data set 
the difference between various fits will be small in this 
range but large in the extrapolated region at the lower tail 
of the distribution where the probability of failure 
necessary for design is located. For simple stress states 
the distribution of flaws may be obtained without assuming 

Q 

any functional form. Matthews* theory was extended for 
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biaxial configurations by Evans. ^ However for multiaxial 

states of stress, a flaw density distribution must be 

extrapolated for untested regions. ** 

To further apply Veibull theory to polyaxial states of 

stress, the Principal of Independent Action (PIA) was 
12 

developed. The PIA hypothesis states that the probability 
of survival is equal to the product of the survival 
probabilities for each of the three principal stresses 
individually. The probability of failure is given as: 


P^ = 1 - exp 



1.2 


where er^, a 2 and are the principal stresses. This method 

will yield unconservative estimates of stress because all 

nonzero principal stresses contribute to the stress normal 

13 

to the plane of the flaw. 

The Batdorf model, ^ accounts for this by including a 
stress space integral inside of Veibull’ s volume/surf ace 
integral. Batdorf ’s assumption that these flaws are cracks 
will combine Veibull’ s statistical model with fracture 
mechanics theory. Material failure is now associated with 
crack growth. Since both normal and shear stresses 
contribute to crack growth, an appropriate fracture 
criterion is needed to account for combined loading. The 
presence of shear reduces the normal stress needed to 
produce fracture. 
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This occurs when the stress on the crack reaches a 

critical value which is a function of the fracture criterion 

and the crack configuration. Two of the mixed-mode fracture 

criteria which have been used in combination with the 

Batdorf model are the maximum tensile strength and the 

strain -energy release rate, for volume cracks and for 

17 

surface cracks. 

Arbitrarily stressed components may be analyzed by 

dividing them into small volumes of material whose stress 

state is assumed constant. Stress analysis using the finite 

element method is compatible with this model to determine 

18 

the failure probability. 

C. Brittle Fracture In Compression 

In general, the models mentioned thus far do not 

include the effect of compressive normal stresses on the 

crack. However, Griffith introduced the idea that brittle 

materials fracture in compression. In a material containing 

pre-existing cracks, the unequal principal compressive 

forces generate shear stresses which act against frictional 

forces producing tensile stresses near the crack tip.^^ The 

crack branches nearly parallel to the direction of maximum 

compression. This secondary crack will grow until the 

tension at the crack tip has dropped to the applied 
20 

compression. This effect is illustrated in Fig. 3. 

Since the existing cracks are microscopic, a single 
crack does not produce total failure as it almost always 
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SINGLE CRACK 



SECONDARY 

CRACK 

ARREST 


Fig. 1.3 Crack growth under compressive loading. 
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does in tension. Total fracture occurs when several of 
these cracks extend and join together creating a shear 
fault. Compressive failure is dependent on the total flaw 
distribution not on the weakest flaw alone. Using the 
weakest link concept, initiation of fracture of a single 
crack has been predicted, for specific stress states 
however, crack interaction was not considered. 

D. Thesis Outline 

The purpose of this study is to model the behavior of 
brittle materials under arbitrary loading conditions. This 
program uses data from a finite element analysis to 
determine the probability of failure originating from volume 
and/or surface flaws. Due to the nature of the fracture of 
brittle materials under compression, the component is 
modeled as a series system in order to establish bounds on 
the probability of failure. When the material is 
represented in this manner crack -crack interaction is taken 
into account. 

In Chapter II, the Batdorf model for failure prediction 
in tension is presented. The material is assumed to contain 
a distribution of uniformly distributed and randomly 
oriented cracks. The failure strength of a material with 
pre-existing cracks is determined from experimental data 
using the 4-point bend test. 

In Chapter III, the modified Batdorf model for 
compressive loading is developed. The reduction in shear 
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due to the effect of the compressive stress on the crack 
face is included in the analysis for volume and surface 
f laws . 

The calculation of the bounds on the probability of 
failure is described in Chapter IV. An element is 
equivalent to a component in the system. Determination of 
the probability of failure for an element is shown based on 
the stress output from the finite element analysis. 

In Chapter V, the model is applied to determine the 
probability of failure for contact stress problems. The 
thesis concludes with Chapter VI, where the results and 
conclusions of the study are presented. The needs for 
future research are given. 


CHAPTER II 


FRACTURE PREDICTION 


The failure of brittle materials has been attributed to 
the presence of flaws. The material fails when the strength 
of the weakest flaw is exceeded. These flaws were assumed 
to be cracks whose strength was dependent on their size and 
orientation. Batdorf assumed that these cracks were 
uniformly distributed and randomly oriented as shown in 
Fig. 2.1. The material will fail when the effective stress 
on a crack reaches a critical value < t qx characteristic of 
that particular crack. To determine the effective stress 
<r e , the shape and the fracture criterion must be assumed. 

Two of the mixed -mode fracture criteria which have been 
used with the Batdorf model are the maximum tensile strength 
and the strain-energy release rate. The strain-energy 
release rate was selected for this analysis because of its 
greater degree of shear -sensitivity . The effective stress 

1 fi 

applied on a crack is given as: 


12 



L3 



NORMAL AND SHEAR STRESSES 
ON AN ISOLATED MICROCRACK 


ibution in brittle materials. 
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a 


e 



+ r 


2.1 


for a Griffith crack, where <r n is the macroscopic tensile 
stress normal to the plane of the crack and r is the shear 
parallel to it. A crack will fracture when the effective 
stress is greater than the critical stress. In the <r -t 
plane eq. 2.1 is a boundary outside of which a crack will 
initiate fracture. This boundary is referred to as the 
fracture envelope. 

A. Surface Flaw Analysis 

Vhen a crack is lying on the surface of a material, it 
is assumed that stresses present only at the surface 
contribute to its growth. As a result the crack is 
subjected to plane stress conditions. The two principal 
stresses on the surface are and ^ 3 a » because they will 
represent the maximum and minimum principal stresses on the 
plane. The stress perpendicular to the surface will not 
effect the growth of surface cracks, however it is included 
in the volume flaw analysis. The subscript a implies a 
surface quantity. 

To determine the probability of failure, the surface is 

divided into elements. Vithin each element the stress state 

is assumed constant. The probability of failure for surface 

24 

cracks may be expressed as: 
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P f a “ P la P 2a 


2.2 


where P^ is the probability of existence in the element of 
surface area of a crack having a critical stress in the 
range <r cr to <r cr + d<7- cr and is the probability that a 
crack will be oriented so that the effective stress is 
greater than or equal to the critical stress. P^ a has the 
form 


P la “ iA 


d V'cr) 

— 5 — — d a 

dt7 cr cr 


2.3 


where A A is the area of the element and N («r ) is the 

surface crack density function which is defined as the 
number of cracks per unit area having a critical stress 
greater than or equal to ^ cr * The crack density function is 
material dependent and is assumed to be 



k Ba'cr‘ 


2.4 


where kg a and m a are determined experimentally. 
P 2a Si ven as: 


P 


2a 



2.5 


where u is the radian measure of the angular range in the 
quadrant of stress space within which the effective stress 
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is greater than or equal to a . The total angular range in 


cr 

this region is 5 . When a > cr over the entire range, 

jL 6cL " cr 


then u = 


5 - 


If a < a everywhere, then u> = 0. 
63> cr 


The 


probability that a crack will be oriented so that cr > a , 
is o < P 2a < 1. 

The overall probability of failure is: 


fa 


= 1 - exp 


-J“ I 


r la dN a (<r cr ) _ 




u d a 


cr 


cr 


2.6 


for the entire area A. The limits on the integration of a 


cr 


are a consequence of the assumed crack density function and 
the stress state. When <T Qr < 0 the quantity P^ a as given in 
eq. 2.3 is equal to zero. ? 2 a is zero if <r cr > because 
the effective stress is never greater than when all of 
the principal stresses are tensile. 

When the crack density function is known, the only 
quantity needed to calculate the probability in eq. 2.6 is 
a function of the stress state and the fracture 
criterion. For a crack on the surface of the material, the 
fracture criterion given in eq. 2.1 is stated as; 


^ea 


~2 5 

a + r 
na a 


2.7 


where is the effective stress acting on a surface crack 

and a and r are the normal and shear stresses on a crack, 
na a 
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The Magnitude of the traction vector |^ a |> is: 

I'al 2 = 'ia cos2 ‘ 9 + 'Sa 51 "^ 2 ' 8 

where /? is the angle between the 0 \^ a -axis and the crack 
normal . 

The normal and shear stresses are respectively given 

as : 


ff na “ 'la" 5 P + '3a sin ' 


2.9 


r 


a 




2.10 


Substituting eqs. 2.9 and 2.10 into eq. 2.7 gives an 
expression for the effective stress, 



( <r la"'3a) cos2 ^ 


2.11 


<r _ > a over the range 0 < /? < , then 

Cx 6a Cr 


Pr..,. ~ COS 
^cr 


-i 


cr 

"2 


*3a 
"2 


^la ” a Zz. 


1/2 


cr 


2.12 


Since u = 0 
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* = 1 'cr^3a 

~w = — cos 

X 

Fig. 2.2 shows the angle u as a representation of the arc of 
Mohr’s circle outside the fracture envelope as given in 
eq. 2.7. 

B. Volume Flaw Analysis 

Vithin an element of volume whose stress state is 
assumed constant, the probability of failure may be 
written 


L/Z 


cr 


- a, 


- a, 


3a 

~2 

3a 


<r 3a <<r cr- <r la 



“ P lv P 2v 


2.14 


where P^ v is the probability of existence in the element of 
volume of a crack having a critical stress in the range a 

CT 

to a QT + d<j- cr and P 2v is the probability that a crack will 
be oriented so that the effective stress is greater than or 
equal to the critical stress. The subscript v will be used 
to define volume quantities. P^ y has the form 


P 


lv 




da 


CT 


CT 


2.15 


where AV is the volume of the element and N (a ) is the 

v ' cr y 

crack density function which is defined as the number of 



L9 



nted by the arc of Mohr’s 
fracture envelope. 
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cracks per unit volume having a critical stress greater than 
or equal to ^ cr * The crack density function is material 
dependent and is assumed to be 


N v = k Bv"cr 


2.16 


where kg v and m v are determined experimentally. 
P 2y is given as: 




2.17 


where 0 is the solid angle on a unit sphere containing all 

of the orientations for which a , > a . The overall 

ev — cr 

probability of failure is: 


fv 


= 1 - exp 


-J " J 


r lv dN v (<r cr ) 


W 


H da 


cr 


cr 


2.18 


for the volume V. The limits on the integration of <r cr are 
similar to those used in eq. 2.6. The only quantity needed 
to calculate the probability of failure is II, a function of 
the stress state and the fracture criterion. For a crack 
inside the volume of material, the fracture criterion given 
in eq. 2.1 is stated as: 


ev 


-f 


nv 


+ r 


2.19 
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where <r ev is the effective stress acting on a volume crack 
and a and r y are the normal and shear stresses on the 
crack . 

A normal to the plane of the crack is defined as shown 
in Fig. 2.3, in principal stress space where 
<r lv < <7 - 2v < ^ 3y . a is the angle between the normal and the 
<r 2v -axis and is the angle between the normal and the 
<r lv -axis in a plane perpendicular to the intermediate 
principal stress. The direction cosines are 


1 = sin a cos 
m = cos a 
n = sin a sin 0 


The magnitude of the traction on the plane is; 


, ,2 2 -| 2 2 2 2 2 

\ a y\ = 'lv 1 + a 2v m + *3v n 


2.20 


2.21 


the normal and shear are 


a nv m -iv 1 ' 




2.22 



Substituting eq. 2.20 into eq. 2.21 

i >2 2 / 2 2\ . 2 2 fl / 2 

k v l = * 3v + (^lv-^v) 8111 cos £ + (^2v 



2 

icos a 


2.23 


2.24 
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and similarly for the normal stress in eq. 2.22 

V = *3v + ('iv-'Sv)*™ 2 * cos2 0 + (^2v-^3v) cos2tt 2,: 

Originally, 1 ^ P 2v was determined in the following 
manner. The fracture criterion in eq. 2.19 may be stated 
using the traction vector 


ev 



then eq. 2.24 is of the form 

a£ 2 + b£ + c = 0 


2.26 


o 

where 6 - cos /3, then 

a = 0 

b = ('lv-'Sv) Sin a 

, 2 2x 2 2 2 

c ” (*2v _<r 3v) cos a + 'Sv " ^cr 

Using the quadratic formula, 6 is determined. Let 
= cos -1 V b then 

H = 1 f sina da 

* Jq 


2.27 
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If 6 < 0, ~0 = ^ which indicates that the effective stress is 
greater than or equal to a everywhere for the angle a. 

V X 

Also if 6 > 1, 0 = 0 or < o cr everywhere for the angle 

a. 

The calculation of ? 2 y i n e Q* 2.27 is simplified if two 
planes tangent to the unit sphere at A and B are defined as 
shown in Fig. 2.4. Their normals are 

n^ = sin a i + cos a j 0=0 

2.28 

ng = cos a j + sin ok 0 = x/2 

where i, j and k are the unit vectors in the ^ v » ant * 

Cg v directions, respectively. The traction on A and B 

, ,2 2 / 2 2-v 2 

I 'A I * ff lv ‘ Wv’^v) 005 “ 

2.29 

, ,2 2 / 2 2x 2 

Kb I = "3V + (^v-'sv) 008 a 
then the traction as given in eq. 2.24 is written: 

k v l 2 = K b I 2 + (KaI 2 - KbI 2 ) 00 ® 2 ' 9 2 - 30 


The normal stresses on A and B are 
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"nA * "lv ‘ ( <7 1 v- <7 2v) cos2 “ 

2.31 

^nB = ff 3v + (^v'^v) 0082 " 
and the normal stress on any plane is: 

'n = ff nB * ( ff nk~ ff ni) cos ^P 232 

P 2v may be determined by substituting eq. 2.30 into 
eq. 2.26 

= I^Bl 2 + ( I 1 2 " I <r B| 2 ) cos2 > 3 cr 2 * 33 

then 



cos 



11/2 


2.34 


where 0 CT is a function of a CT , |<r A |, |<r fi | and implicitly a. 

Eq. 2.34 is similar to eq. 2.12 for surface cracks. As was 

the case for surface cracks, cr > a over o < 3 < 0 , or 

ev — cr — — ^ C r 

u - P CT ‘ E<1* 2.27 may be written: 


n 


t/2 _ 

sin a da 


= a 

J 0 


2.35 


where 


27 


u = 0 


U! = 


I cos- 1 


cr 


B' 


Bi 


1/2 


•A I ^ ff cr 


'Bl < a cr < \*k\ 


2.36 


u = 1 


'cr * I 'B 


where |<7^| and |<7g| are given by eq. 2.29 as a function of 
a. The region inside Mohr’s circle of stress which lies 
outside the fracture envelope represents orientations within 
which a crack must lie to initiate fracture as shown in 
Fig 2.5. An example of the solid angle 0 on a unit sphere 
is shown in Fig. 2.6. 


C. Determination Of Material Parameters 

The crack density functions, N a and N v , contain two 
constants which must be determined experimentally. A plot 
of the probability of failure versus the fracture stress is 
made from the data and is fitted to: 



where is the fracture stress, k and m are determined 
using the least squares method. 

Given a fracture theory, an expression for the 
probability of failure can be found which is a function of 
the specimen geometry and loading conditions. This 
expression is used in combination with eq. 2.37 to evaluate 
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the material parameters. Ideally, a uniform tension test 
would be the easiest to formulate. However, this test is 
not popular because the ceramic tensile specimens are costly 
and their load train is not easily aligned. The four -point 
bend test has become the preferred test because it can be 
controlled and the stress field is uniaxial but not uniform. 
The test configuration is shown in Fig. 2.7 for a beam with 
a circular cross-section. 

The stress distribution in the beam as a function of 
the maximum tensile stress is given by 

x r sin 6 0 < x < (L q -L^)/2 

. Q (L o -L i )/2 < x < (L Q +L.)/2 2.38 

(L q -x) r sin 6 (L Q +L i )/2 < x < L q 

where E. is the radius of the beam, L Q is the distance 

between the outer loads, is the distance between the 

inner loads and a is the magnitude of the maximum tensile 

max 

stress in the beam, a - a „ between the inner loads and at 

max 

0 = Since the magnitude of the tensile stresses and the 
compressive stresses are equal, the beam will most likely 
fail in tension. 
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Fig. 2.7 A beam of circular cross-section with four point 
loads . 
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1. Surface Flaws 

For a uniaxial state of stress with <r^ a = <r, and 
0 g a = 0 the probability that a crack will be oriented so the 
effective stress is greater than or equal to the critical 
stress is: 


U) 


2 

r 


cos 


-1 



substituting eq. 2.39 into eq. 2.6 


2.39 


P fa= 1 


exp 


-i 1 1 

A 0 


cos 


-1 


cr 


dN aKr) 

— | — — da 
d<r cr cr 


dA 


2.40 


If it assumed that when a = 0, then N = 0, eq. 2.40 is 

cr ct 

integrated by parts and substituting eq. 2.4 into eq. 2.40 


fa 


= 1 - exp 


■J.J 


A J 0 ~ 


m a 

2k fi (r a 

° r 2-j 1/2 da cr dA 
cr-' 


2.41 


Let 


V = 



then 
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P fa ■ 1 - e *P | -J 


k Ba <T 


Z dA 


2.42 


where 


z -J * 

0 


K- 1 

T3 - 


{I -,}’ 1 / 2 d„ 


2.43 


Integrating eq. 2.42 over the tensile portion of the 


beam 


P fa = 1 - «p 


in 


7T L n’t a 

o Zk Ba o- 


•I J 


Rdxd(? 


0 0 


2.44 


and substituting eq. 2.38 into eq. 2.44 gives an expression 
for the probability of failure in terms of the specimen 
geometry and the loading conditions 


P fa = 1 ' ex P 


ZKk Ba(V"a L i) "i 


m 


a t 


max 


*( m a + l) 


r 

J sin 


e d 9 


2.45 


The constants kg a and m a are determined by equating the 
exponents in eq. 2.45 and eq. 2.37 


m = m 
a 


2.46 
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5r(m +l)k 

Jr _ S 

Ba rir m 

ZE(L +m L. ) sin a 9 d 6 
o a i Jq 

where a is the fracture stress, <r, . 
max 7 I 


2.47 


2. Volume Flaws 

For a uniaxial state of stress with a ^ = er, and 

a 2 V = tr^y = 0 the probability that a crack will be oriented 
so the effective stress is greater than or equal to the 
critical stress can be computed directly^ 


n = i - 


cr 


2.48 


substituting eq. 2.48 into eq. 2.10 


p fv = 1 - exp 


V 

-n 

V 0 


cr 

a 


dN (a ) 

d 'cr 

cr 


2.49 


substituting eq. 2.16 into eq. 2.49 and integrating over cr 


cr 


P fv = 1 - exp 


-J 


k Bv _“v HV 


2.50 


Integrating eq. 2.50 over the tensile portion of the 


beam 
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R L, 


P fv = 1 


exp 


•j i j 

0 0 0 


k n m 

7 — \ r a v r dx dr dO 

v m v +1 ; 


2.51 


and substituting eq. 2.38 into eq. 2.51 gives an expression 
for the probability of failure in terms of the specimen 
geometry and the loading conditions 


P f v = 1 - exp 



9 

(L +m L.) Rr a v 
v o v i y max 

(m v +l) 2 (m v +2) 



6 d0 


2.52 


The constants kg y and m y are determined by equating the 
exponents in eq. 2.52 and eq. 2.37 


m y = m 


2.53 


(* v +l) 2 («v+ 2 ) k 


L Bv 


9 r x m v 
(L +m L.) IT sin v Q d0 
K o \ \> J 0 


2.54 


where <j is the fracture stress, er , . 
max x 


D. Numerical Integration 

Vhen the stress state and the crack density function 
are known, the probability of failure may be calculated. 
However, numerical integration is necessary to compute the 
failure probabilities as given in eqs . 2.6 and 2.18. To 
evaluate these integrals the Gaussian quadrature method is 
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used. The integration of an arbitrary function over a 
finite interval may be approximated by 



where N is the number of sampling points, x^ and w^ are the 
location and weight of the ith sampling point, respectively. 
These values x^ and w^ are tabulated for different numbers 
of sampling points. 

The probability of failure for surface cracks is more 
efficiently computed^ if 


S 


a 



then substituting eq. 2.4 into 2.6 gives: 


2.56 



m r 

1 m -1 

P fa = 1 - exp 

‘ Am a k Ba (r la J 

u S * dS 0 
o a a 


where ~w is given by eq. 2.13 or as given by eq. 2.36 where 
|<r^[ = <r^ and |<7g| = Applying eq 2.55 to the integral 

in eq. 2.57 gives: 


-Am_k D _ a 


m_ 


N 


a Ba la 


m_ -1 


I "i w i 

L X cl • X 

i=l 1 


p fa * 1 - ex P 


2.58 
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where 




+1 


7 


u- is evaluated at £ and x„ and w. are tabulated 
i *a^ a^ i 

according to the number of sampling points used. 

A similar process is utilized to evaluate the 

probability of failure for volume cracks as given in 

eq. 2.18. If 


S 


v 



2.59 


then substituting eqs. 2.59 and 2.16 into eq. 2.18 gives: 



m r^/ 2 r 

1 m -1 

p fv = 1 - ex P 

- V »v k Bv'lv J Q J 

w sin a S dS da 

0 v v 


2.60 


where u is given by eq. 2.36. Applying eq 2.55 to the 
integral in eq. 2.60 gives: 


p fv = 1 " exp 


m 


" N 


IT 


P r ®y “ 1 

) ) u-- sin (. £ w.w. 

L L 13 s j s v. 1 3 

i=l j=l 1 

2.61 


where 
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V 1 


ir(x.+l) 

C j 1 ~ 


w.. is evaluated at £ and £. . x , w. , x. and w. are 

J ^ J ^ J J 

tabulated as a function of M or N. 

It is of interest to note that for both volume and 
surface analysis the probability of failure is dependent on 
a common function, w. Therefore in the subsequent analysis 
of compressive stress states, it is not necessary to 
consider both volume and surface probabilities but to 
formulate ~w so that it may be used in both eq. 2.58 for 
surface probability and in eq. 2.61 for the probability of 
failure for an element of volume. 



CHAPTER III 


FRACTURE INITIATION UNDER COMPRESSIVE LOADING 


Brittle materials with pre-existing cracks may fracture 

21 

when loaded in compression. The unequal principal 
compressive stresses generate shear stresses which act 
against frictional forces, initiating local crack growth. 
Since the existing cracks are microscopic, a single crack 
does not produce total failure as it almost always does in 
tension. Total failure occurs when several of these cracks 
extend and join together. Before the failure of an entire 
component can be analyzed, the fracture of a single crack 
must be considered. 

Vhen the stress normal to the crack plane is 

compressive, the shear and the friction due to the normal 

23 

stress will act against one another 

r e = |r| + ti<r n 3.1 

where r is the effective shear stress on the crack and /jl is 
e 

the internal friction coefficient of the material. r e is 
defined as the shear stress necessary to initiate fracture 
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when the normal stress is less than zero. 

The criterion given in eq. 3.1 is combined with the 
strain-energy release rate criterion in eq. 2.1. must be 

equal to the value of the effective stress as given in 
eq. 2.1 when = 0, in order to preserve continuity in 
shear. Then the fracture criterion for < 0 is stated: 

= T e ~ M + 3 ’ 2 

The fracture envelope in shown in Fig. 3.1. Failure 
probabilities are based on the combination of the two 
criteria as were given in eqs. 2.1 and 3.2. 

To accomodate the new fracture criterion, two 
modifications are made to the Batdorf model as presented in 
Chapter II. First, the maximum effective stress ^ emax > is 
not always equal to the maximum principal stress as it did 
when all of the principal stresses were tensile, resulting 
in a change in the limits on the integration of or in 
eqs. 2.6 and 2.18. Second P 2 the probability that the 
effective stress is greater than or equal to the critical 
stress, must be reformulated. 

The subscript notation used in this and subsequent 
chapters will not include a reference to a surface or volume 
flaw analysis. For example, the quantity N refers to both 
N the surface crack density function and N the volume 

ct ’ 

crack density function. 




Fig. 3.1 Fracture envelope. 
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A. Maximum Effective Stress 

In order to calculate the effective stress exceeding 
the critical stress P 2 , the effective stress for all crack 
orientations must be considered. Vhen the critical stress 
is greater than the maximum effective stress, P 2 = 0. In 
terms of the fracture envelope shown in Fig. 3.1, a = 

6 cIHcLX 

when Mohr’s circle of stress is tangent to the fracture 

envelope. The orientation of the normal to the plane on 

which <j acts must lie in the 13 -plane, 
emax 

The maximum effective stress is characterized 
differently depending upon the nature of the stress on the 
element . 

1. Compression - Mohr’s circle is tangent to the linear 
portion of the fracture envelope 

2. Tension - The maximum effective stress is equal to 
the maximum principal stress . 

3. Tension and compression combined - A transition 
region where ^ emax is equal to the shear stress when 
only shear is present. 

The limitations of each of these three methods will be 
discussed in order to develope a general algorithm for the 
determination of the maximum effective stress. 
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1. Compression 

Mohr’s circle is tangent to the linear portion of the 
fracture envelope when the effective stress given by eq. 3.3 
is a maximum. The normal and shear stresses in the 13 -plane 
are : 


°l +0 3 a l~ a 3 


a n = 2 + 2 cos 




r | = — ^ — sin 2 /3 


3.3 


where 0 is the angle between the crack normal and the 
o^-axis. The effective stress in the 13-plane is given by 
substituting eq. 3.3 into 3.2 


a l~ a 3 


sin 2 0 + 


A 


a l +a 3 

— t r~ 


— cos 2/? 


3.4 


To determine where a is located, one must find the 

emax ’ 

orientation where the slope 

da 

— - = (a.-ao) cos 20 - (i (a.-a?) sin 20 
80 ° . 

is equal to zero. Solving for 0 = /? v , at a = <r one 

niajc 0 emax 

obtains 


tan ^max 


1 

A 


3.5 
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which is a function of the internal friction coefficient. 
Mohr’s circle is tangent to the fracture envelope at the 
point where 0 = /3 „ as shown in Fig. 3.2. For the limiting 
case ji = 0, ^ emax occurs at P max = which is the location 

of the maximum shear stress. 

A general expression for ^ emax can be found from 
eq. 3.5, given by: 


sin 2 0 


max | k 

Jl + ft 


cos 20 = ^ 

■max | 7 ) 

Jl + / 


substituting eq 3.6 into eq. 3.4 gives: 


3.6 


a 


emax 


J 1 


a l' a 3 

— 2 ~ 


+ [I 


2 


3.7 


Since a ^ > a the first term in eq 3.7 is always greater 
than or equal to zero. However <r < 0 for all 0 , then 
a* + a o < 0. If <r <0, the frictional force along the 
crack is greater than the shear and all cracks under these 
loading conditions are locked. For example, under 
hydrostatic pressure cr^ = < 0, no shear is present. This 

phenomenon also takes place near the orientation of 
compressive principal stresses. 

For any stress state whose maximum effective stress is 
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less than or equal to zero, cracks will always lock. The 
probability of failure for these regions is always equal to 
zero, therefore they may be eliminated from the overall 
analysis. This is done by considering the ratio of minimum 
to maximum principal stress: 



when Og < 0 ^ < 0, R > 1 because the magnitude of is 
greater than or equal to the magnitude of When 

£r =0, R = IU then from eq. 3.7 

(1+p 2 ) 1 / 2 + it 

*° = y/s . „ 3 - 8 

Fig. 3.3 shows the crack locking region in the <r -t plane 
with Mohr’s circle at R = R Q . 

Rq as a function of the internal friction coefficient 
is shown in Fig. 3.4. If 1 < R < Rq the maximum effective 
stress is less than zero and all cracks will lock. The 
elements whose principal stress ratio falls in this range 
are eliminated from the failure analysis. 
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Fig. 3.3 Mohr’s circle vhen the maximum effective stress is 
equal to zero. 
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2. Tension 

When the normal stress is greater than or equal to zero 
everywhere or 0 < < <r. , 0 - is determined using the 

strain-energy release rate fracture criteria. As was 
discussed in Chapter II, the maximum effective stress is 
equal to the maximum principal stress. Mohr’s circle is 
tangent to the fracture envelope at a n = a ^ as shown in 
Fig. 3.5. 

3. Tension and compression combined 

The maximum effective stress has been determined if the 

normal stress is either tensile or compressive for any crack 

orientation. When <r^ < 0 < tr^, both of the fracture 

criteria must be considered to find «r „ . The value of 

emax 

^emajc is S°verned by: 

a. Compressive criterion. 

b. A transition region between tensile and compressive 

criteria. 

c. Tensile criterion. 

depending on the ratio of minimum to maximum principal 
stress . 

a. Compressive criteria. 

The maximum effective stress will be located at 

/? = /? v as long as the normal stress at that orientation is 

compressive. The normal stress at the orientation of 

maximum effective stress, = 3_ 

max 


is given by eq. 3.3: 
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°n = 


*1**3 


<r l~ (r 3 


+ — 5 — cos 20 ] 


max 


3.9 


To determine when the normal stress changes from compressive 
to tensile, let a = 0 and substitute eq. 3.6 into 3.9 which 
gives : 




(1 +/i 2 ) 1/2 + p _ 

(i +/ , 2 )V2 . ^ ■ ‘"o 


3.10 


where Rq is defined as the principal stress ratio when 
^emax = If R £ -Rq the maximum effective stress is given 
by eq. 3.7 because < 0 at /? m , or Mohr’s circle is 

tangent to the linear portion of the fracture envelope as 
shown in Fig. 3.6. 


b. The transition between tensile and compressive criteria 

If the normal stress at 0 = /? is tensile the maximum 

max 

effective stress may occur at the angle which separates the 

tensile and compressive fields or <r n = 0. Let 0 = 0q at 

<r = 0 from eq. 3.6 
n 


2/? n = cos 


-1 


*1**3 


a l~ a Z 


3.11 


then the shear is: 


emax 


" T 0 " 


'1 ~°3 

■V 5 sin 2/3 q 


3.12 
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where Tq is the magnitude of the shear stress at <r n = 0. 

Fig. 3.7 shows an example of a state of stress where 

a = r A . The maximum effective stress will be equal to 
emax 0 

Tq as long as t q > a^. r o = if ~ ~ a 3> f* 1611 
<r = r A over the range: -R n < R < -1. 

c. Tensile criterion. 

Vhen the tensile stress dominates, the maximum 
effective stress is equal to as it does when a • is always 
tensile. For | | < <r^, ^ emax = ^ as shown in Fig. 3.8. A 

summary of the maximum effective stress for different stress 
states is given in Table 3.1. 
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Table 3.1 The maximum effective stress for ratios of 
principal stress. 
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B. Fracture Prediction 

The probability of failure for an element whose stress 
state is assumed constant is given by eq. 2.57 as: 

P fa = 1 ' exp 

for surface elements and by eq. 2.60 as: 



m rx/2 r 

1 _ m -1 

1 - exp 

-Vm kn cr 

v Bv emax J Q J 

u sin a S v dS„ da 

q cr cr 


3.14 


-Am k n cr 

a Ba emax 


m a f 1 - 

U) 

ax J Q 


m a _1 

S a dS 
cr cr 


3.13 


for volume elements, where <r is replaced by a , then 

1 J wmax ’ 


cr 


cr a. 


emax 


and u is the probability that a crack will be oriented so 
that the effective stress is greater than or equal to the 
critical stress formulated using the appropriate fracture 
criterion. H) will be a function of two angles 0^ and f}^ and 
is expressed as: 


w = 


3.15 


As is the case for the maximum effective stress, u is 
dependent on the nature of the stress on an element and is 
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different for the normal stress which is compressive, 
tensile or compressive and tensile combined. 


1. Compressive loading. 

H is reformulated using the fracture criterion for 
<r n < 0, given in eq. 3.2. Substituting 


T 


2 



into eq. 3.2 gives: 


a 


e 



+ »<r T 


3.16 


Using the notation given in Table 3.2 for surface or volume 
analysis, eq. 3.16 becomes: 


S e * 


l 5 + 


fj. s 


3.17 


which may be written: 

(S e - /,S) 2 = T - S 2 


or 


S 2 - 2/iS S + (lV)S 2 - T = 0 

V C 


3.18 
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Table 3.2 Notation for volume and surface flaw 
analysis . 



Critical Stress 


Critical Stress (Squared) 


Effective Stress 


Normal Stress at A 


Normal Stress at B 


Traction at A 


Traction at B 
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The traction and the normal stress are: 

T = T fi + (T a -T b ) cos 2 /? 3.19 

S = S fi + (S A -S fi ) cos 2 /? 3.20 

Substituting eqs. 3.19 and 3.20 into eq. 3.18 gives: 

S e ■ 2/iS e S B - 2/zS e (S A -S B )cos 0 

+ (1 +[i )S B + 2(1+// ) S b (S a ~S b ) cos 0 + (l+/x )(S A -S fi ) cos (3 
- T fi - (T A -T B )cos 2 y5 = 0 3.21 

At S = S , 0 = 0 and eq. 3.21 may be rewritten: 

C vX v X 

o 

a 6 + b6 + c = 0 

o 

where 6 = cos /? cr , then 
a = (1+// 2 ) (S A -S fi ) 2 

b = -2//S cr (S A -S B ) + 2(1+// 2 )S b (S a ~S b ) - (T a -T b ) 
c = T cr - 2//S cr S B + (1+// 2 )S b - T b 

then using the quadratic formula: 


6 = 


-b ± 



- 4ac 
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where 6 has two roots, and 6 n . S„> S„ over the 

7 l 2 e _ cr 

interval < /? < Pq2' ^ Si yen as: 


2 ^C2'^C1^ 

u = 

X 


3.22 


where 


w = 0 
/?C1 = o 

/?C2 = */2 


b 2 - 4ac < 0 

6 1 > 1 

^2 * 0 


then 0 < &; < 1. The region inside Mohr’s circle and outside 
of the fracture envelope represents orientations within 
which cracks must lie in order to initiate fracture as shown 
in Fig. 3.9. An example of the solid angle H on a unit 
sphere for volume distributed cracks is shown in Fig. 3.10. 
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2. Tensile loading 

The probability of failure for surface and volume 
elements under tensile loading is discussed in Chapter II. 
To use the general formula for Hi given by eq. 3.15, let 
0^ = 0 and 02 = P i v h ere 


0r^ - 0 

/?p = cos 
= jt/2 



B 

B 


T A * T cr 


T B < T cr < T A 


3.23 


T cr i h 


T. , T n and T are defined in Table 3.2 and /? T = /? as 
A’ d cr 1 cr 

given by eq. 2.12 for surface elements and eq. 2.34 for 
volume elements. 


3. Combined tensile and compressive loading 

In order to compute ~u for elements whose loading is 
both tensile and compressive, the two fracture criteria are 
needed. The tensile criterion as given by eq. 2.26 is 
expressed using the notation in Table 3.2 as: 

S e = T 1 / 2 = [T fi + (T a -T b )cos 2 ^] 1 / 2 3.23 

and the compressive criterion is given by eq. 3.17. The 
range over which the effective stress is greater than or 
equal to the critical stress is 0^ < 0 < @ 2 > v h ere 0-^ P% 
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are equal to the appropriate combination of 0^ t 0 q^ or 0q%' 
To determine that combination each criterion must be 
considered separately. The effective stress, S e , as a 
function of 0 is shown in Fig. 3.11a for the tensile region 
and Fig. 3.11b for the compressive region. 

0 ^ and 02 are dependent on the maximum effective stress 
for either criteria and whether the normal stress is tensile 
or compressive at that point. For any orientation ( a,0 ) the 
maximum effective stress is located in the 13-plane, 
a = tt/2. Vhen a f ;r/2, the maximum effective stress is not 
so easily defined. For the tensile criterion the maximum 
effective stress is equal to T^ 1 / 2 or Tg 1 / 2 , whichever is 
larger. 

For the compressive criterion, the maximum effective 
stress is determined using Fig. 3.12. Fig. 3.12 shows 
Mohr’s circle relative to a line whose slope is the same as 
the linear portion of the fracture envelope. By definition, 
Mohr’s circle for the 13 -plane has its center at point 0, or 
the average of the minimum and the maximum principal 
stresses. As a is constant from A to B, (shown in Fig. 2.4) 
the normal and shear stress (<r n ,r) comprise an arc of a 
circle which is concentric with Mohr’s circle for the 
13 -plane. The stress (<r n ,r) at the orientation of maximum 
effective stress will be located on the line UP as shown in 
Fig. 3.12 or the point on the arc closest to the line UP. 
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/S (rad) 

a) TENSILE CRITERION 



Fig. 3.11 The effective stress as a function of crack 
orientation for A)tensile and b)compressive 
loading. 
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T 



/ 


Fig. 3.12 Mohr’s circle with the shear as a function of the 
normal stress for constant a. 
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The normal stress along OF as a function of a is: 

a = R cos 2 0 + — n — 3.24 

nm a max z 

where R is the radius of the circle on which the arc AB 
a 

lies and cos 2/? max is given by eq. 3.6. If a nfi < er^ < <r nA , 

the maximum effective stress is located where = <r which 

n nm 

is on the line UF. The shear is: 




sin 20 


max 


3.25 


where sin 20 is given by eq. 3.6. The maximum effective 

H13X 

stress is found by substituting eq. 3.24 and 3.25 into 
eq. 3.2 


a = 
em 


1 + 


"7 K 


a l +a 3 


3.26 


B. = (a^-cr 0 )/2 if a = t/2 then = a for the 13-plane 

as given by eq. 3.7 If <r nB > g , n|B , the maximum effective 
stress occurs when <r n = <r nB and when is a 

maximum when a = <r nA . For "the example shown in Fig. 3.12 
as a approaches zero, <r nB > the maximum effective 

stress occurs at 0 = t/2. 
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The angles 0^ and 0^ are defined differently for: 


- > Tg and <r > 0, The traction is decreasing as 0 
increases and the normal stress is tensile when the 
effective stress is maximum in compression. 

- T^ > Tj and <r < 0, The traction is decreasing as 0 
increases and the normal stress is compressive when the 
effective stress is maximum in compression. 

- T^ < Tg and > 0, The traction is increasing as 0 
increases and the normal stress is tensile when the 
effective stress is maximum in compression. 

- T^ < Tg and < r nm < 0, The traction is increasing as 0 
increases and the normal stress is compressive when the 
effective stress is maximum in compression. 


An example of the first case T^ > Tg and <r > 0, is 
shown in Fig. 3.13. Let S e = S cr , then 0^ and 0 2 are: 


X* 

II 

O 



O 

II 

CS 

t a 1/2 

< s 

- cr 

II 

cs 

T 0 1/2 

< S 
- cr 

@2 ~ &C2 

S cr < 

X 

L 0 



where T = Tq at 0 = 0q. 
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The second case > Tg and < 0, is physically 
impossible. Consider the fracture criterion for S < 0: 


S 


e 


T - S 2 + /zS 


if the traction and the normal stress are both decreasing 
then the effective stress must decrease. Vhen the traction 
is increasing the effective stress may have a maximum when 
S < 0. 


The effective stress as a function of 0 for the third 
case T. < T n and a > 0, is shown in Fig. 3.14. Let S = 
S cr , then 0^ and 0 ^ are: 

h 
h 
@2 
02 


0 

0r£ 

0 

0 

0Q2 


V' S S cr 

T A V2 < S cr < T 0 

S cr i T A V2 
T 0 1/2 < S cr 


1/2 


S cr < T 0 


1/2 
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FRACTURE CRITERIA 

COMBINED 
• • • • TENSILE 



P (rad) 


Fig. 3.14 The effective stress as a function of crack 

orientation, T, < Tr, and a >0. 

A d run 
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The effective stress as a function of 0 for the final 

case T. < T n and a < 0, is shown in Fig. 3.15. Let S 0 = 

A ” d nm e 

S cr , then 0^ and 0^ are: 

0 . = 0 S pm < S_ 

1 em cr 



— 0^ 

t a 1/2 

^ s < s 

cr em 


ii 

H 4 

s < 

cr - 

t a 1/2 


= 0 

S < 

s 



em - 

cr 


II 

to 

S cr < 

S em 


After computing the probability of failure for an 
individual element under different loading conditions, these 
elements are assembled in order to analyze the failure 
probability for an entire component. 
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FRACTURE CRITERIA 


-e— COMBINED 
• • • TENSILE 



rr/2 


P (rad) 


Fig. 3.15 The effective stress as a function of crack 

orientation, T, < T™ and a < 0. 

A d nm 


CHAPTER IV 


SYSTEM RELIABILITY 


The probability of failure for an element of area or 

volume has been evaluated for tensile or compressive loading 

in Chapters II and III. In compression, several cracks 

extend and join together to create a shear fault as shown in 
20 

Fig. 4.1. To account for the multicracking phenomenon, 
the material is modeled as a series system rather than 
independent elements as in the weakest link theory. 
Reliability analysis is used to correlate the elements and 
establish bounds on the probability of failure. 

A. Finite Element Analysis 

This analysis presented to this point has been 
concerned with evaluating the probability of failure within 
an element of area or volume where the stress state is 
assumed constant. To evaluate these stresses the finite 
element method is used. An element of volume or surface is 
not the same as an element in "finite element", because its 
stress distribution may not be assumed constant. The finite 
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Fig. 4.1 The joining of several cracks to create a shear 
fault under compressive loading. 
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element must be divided into sub -elements so that the 
constant stress assumption is valid. The number of 
subdivisions is governed by the number of locations at which 
the stress is output from the finite element program. 

Along with the stress at each one of these locations, 

the area or volume of that sub -element is determined. 

Finite element analysis makes use of isoparametric elements 

within which the displacements are formulated using 

27 

interpolation functions. These interpolation functions 
enable an element of arbitrary shape in the global 
coordinate system to be mapped into a natural coordinate 
system over which the calculations are carried out. The 
global coordinate system (x,y,z) is representative of the 
physical system. The natural coordinate system (r,s,t) is 
constructed so that -1 < r < 1 , -1 < s < 1 and -1 < t < 1, 
The volume of an element is determined using the Jacobian 
matrix 


[J] 


dx dv dz 

WWW 

dx dv dz 

WWW 

dx dy dz 

WWW 


where for an element with N nodes: 


N 

dx v 

w =2 

i=l 




x i 


4.1 


dx 


4.2 
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then is the interpolation function and x^ is the global 
x-coordinate of the ith node. Each term in the Jacobian 
transformation matrix is computed in the same manner and is 
a function of r, s and t. 

The volume of an element is given by: 
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I I I |m 
-1 -1 -1 


dr ds dt 


4.3 


which may be evaluated numerically using eq. 2.55: 



I J K 


III [J( r i> s r t k )] 

i=l j=l k=l 


w i w j w k 


4.4 


where I, J and K are the number of integration points in the 
r, s and t directions, respectively. The volume of a 
sub-element is: 


V = 


[ J ( r i» s j ,t k )] 


W i W j W k 


4.5 


where w. , w- and w, are the weights associated with each 
i j x 

sampling point. The number of sub -elements will be equal to 
the number of integration points used for the element. This 
technique is appropriate if the stress has been output at 
the integration points. 

The area of a sub -element is determined in the same 


manner : 



where the Jacobian as given by eq. 4.1 is reduced to a 2x2 
matrix. 


B. Reliability Bounds 

Bounds on the probability of failure are determined by 
considering the elements of a component and their 
relationship to one another. Before that relationship can 
be found, the method over which the entire component fails 
is determined. The fracture of a brittle material will 
occur when a crack in any part of the component fractures . 

This weakest link hypothesis is an example of a series 
system which is defined in the following way. The component 
or structure is considered at a fixed point in time, the 
status of that structure (functioning/failed) is dependent 
on the states of its elements. The state of an individual 
elements is expressed in terms of two binary variables, 

A and B: 


A i 


1 if element i is functioning 
0 if element i has failed 


4.7 


where i ranges from 1 to k and k is the number of elements 
in the component. The state of the structure is: 



1 if the structure is functioning 
0 if the structure has failed 


1 - A s 4.8 
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where A g for a series system is a function of all of its 
elements and is given as: 


A 


s 


A^Aj . 


.A^ = min 



Substituting eq. 4.8 into 4.9 gives: 


4.9 


A s “ A l A 2 ,,,A k-l ' A i A 2 • • • A k-l B k 


Repeating this operation yields: 


Bg — + A^2 + + . . . + A^A2* • . Aj^^Bj^ 4.10 


Since the state variables can only take on values of zero or 
one, it follows that: 

A^Ag-.-A^ > max|l - (B^ + Bg + ••• + B^ ),0| 
which when combined with eq. 4.10 leads to: 


B s * B 1 


k 

l 


+ 2 max 

i=2 


i-1 


B i 


I ®i B j »° 


3=1 


4.11 


For any j<i it also follows that 


A ^ An « • • A • ^ A • ™ 1 ” B • 
1 z 1 " J J 


4.12 
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combining eq. 4.10 and 4.12 


b < y b. 

s - L i 

i=l 


V max B . B . 

L i<l 1 J 
i=2 J<1 


4.13 


Bounds on B are given by eqs. 4.11 and 4.13 or 
s 


B l+ 1 

i=2 


max 


i-1 

V 1 B i B j’° 

3=2 


* B s * l B i 
i=l 


l 


max B . B . 4.14 

" 1 J 

i=2 J<1 


Ditlevsen bounds on P* are similar to the bounds on B 

I s 


as given by eq. 4.14 
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I 

i=2 


max 


i-1 

p i- 

j=2 


S P f * I T i 

i=l 


I 


max P. . 

" -i^l X J 

i=2 3<1 


4.15 


where P^ is the probability of failure of an individual 
element, assembled in decreasing order and P. . is the joint 
probability of failure of elements i and j . is equal to 

the probability of failure for the element with the highest 
probability of failure. P.. is formulated for simplicity as 

J 

a function of the safety indices 0 . and 0 - , given by 

J 


’ij = P i P j + J 0 >^j :z ) 


dz 


4.16 


where 
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0>(x,y:p) = exp 

23 r Jl - p 2 

p is the correlation coefficient and the safety index 0^ = 
-<J> -1 (P^) is a normally distributed function. The 
probability of failure as a function of the safety index is 
shown in Fig. 4.2. 

If the correlation coefficient p is equal to zero, the 
elements are not correlated and if p = 1 they are fully 
correlated. When p is equal to one, a k series system is 
modeled as one single element whose probability of failure 
is the average of the k elements. 

In order to minimize the size of the problem, the 
weakest link model is used to combine the sub -elements 
within an element so that the element may be used to 
determine the reliability bounds. 

The probability of survival is: 







CHAPTER V 


APPLICATION TO CONTACT STRESS PROBLEMS 


The preceding theory is used to evaluate the 
probability of failure for an alumina ceramic under two 
different contact stress conditions. First, pressure is 
applied to two cylinders in contact and second, a 
compressive load is applied to a beam with a machined notch. 
The second example is a model of a test conducted at 
NASA Lewis Research Center where a compressive load was 
applied to a sample with a machined notch in order to 
initiate a pre -crack for subsequent studies on fracture 
toughness. For the first example the stress distribution in 
closed form is known. However, in general the 
displacements and stresses are found using the ADINA finite 
element program. 

The material chosen was a 967. alumina ceramic which was 
used in the experiments at NASA. Three -point bend test data 
on a beam of circular cross-section was available in the 
literature and is given in Table 5.1. 
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TABLE 5.1 Fracture stresses and failure probabilities 
for 3 -point bend specimens. 


Fracture Stress 
<r f (MPa) 


Probability of Failure 


P 


f 


378 

417 

421 

430 

448 

453 

455 

457 

461 

470 

472 

475 

479 

493 

495 

497 

502 

528 

532 

540 


0.048 

0.095 

0.143 

0.190 

0.238 

0.286 

0.333 

0.381 

0.429 

0.476 

0.524 

0.571 

0.619 

0.667 

0.714 

0.762 

0.810 

0.857 

0.905 

0.952 


The 4-point bend test is preferred because the shear stress 
between the inner loads is zero. However, the 
length -to -radius ratio for this beam was large enough so the 
maximum shear is an order of magnitude less than the maximum 
normal stress, and is therefore, neglected. The material 
parameters are calculated from this data and are given in 

on 

Table 5.2. The dimensions of the beam are given, k and m 
are calculated using the least squares method and kg a and 
kgv are determined using the technique described in 
Chapter II. 
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Table 5.2 Alumina test data. 

R = 0.0016 m 

L = 0.0254 m 
o 

= 0. (3 -point bend) 

m = 12 . 2 

k = 1.80 x 10" 33 MPa -12 ' 2 

k Ba = 3,99 X 10 " 27 MPa -12 * 2 /m 2 
k fiv = 9.72 x 10" 23 MPa -12 * 2 /m 3 


The last material parameter necessary for this analysis 
is the internal friction coefficient. An exact value for 
this quantity is not known and it must be approximated. The 
internal friction coefficient is defined as being equal to 
the slope of the fracture envelope. When the tensile and 
compressive strengths are known, using simple geometry the 
friction coefficient is calculated. Using the bend strength 
of alumina, 315 MPa (45 ksi) , and the compressive strength, 
2,625 MPa (375 ksi), the internal friction coefficient, /*, 
is determined and is equal to 1.27. 
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A. Two Cylinders In Contact 

Vhen pressure is applied to two cylinders in contact, 

(bearings), the region beneath the load is subjected to high 

compressive stresses. If the displacement in the contact 

region is assumed to be uniform, a closed form solution for 

the stress distribution can be found assuming a 

29 

semi -inf inite region. A schematic representation of the 
applied load on the surface is shown in Fig. 5.1a where Pq 
is the maximum pressure at the center of the contact area 
and a is the half -width of the contact area. The maximum 
pressure is: 



5.1 


where p is the applied load and L is the length of the 
cylinder. The half -width of the contact area is: 


where 






5.2 
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a) The normal and shear load on the boundary 



b) Mesh 


Fig. 5.1 Schematic view of a contact stress distribution 
on a semi-infinite region. 
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where R^, and are the radius, Poisson’s ratio and 
Young’s modulus of the first cylinder and R 2 , ^ are 

respective properties of the second cylinder. 

A closed form solution is available for this loading 
configuration, in order to evaluate the failure probability, 
the stress and the volume/area at discrete points is needed. 
The mesh used to divide the contact region is shown in 
Fig. 5.1b. The stress is evaluated at the center of each 
element. The rectangular area is easily determined, 
however, as the load increases the width of the contact area 
increases. To compensate for this increase, the volume/area 
is calculated as the load increases using eqs. 5.1 and 5.2 



5.3 


where p n is the initial maximum pressure and a. is the 
u i 1 

corresponding half -width of the contact area. The initial 
load along with the bearing dimensions are given in 
Table 5.3. 


C- 3- 
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Table 5.3 Bearing material properties and loading 
conf iguration . 

Dimensions 

R 1 = R 2 = 0.01 m 

L = 0.01 m 

Material properties 

E 1 = E 2 = 280,000 MPa 
= 0.25 

Initial loading condition 
p t = 82.5 N 
a t = 1.875 x 10' 5 m 
Pq = 280 MPa 


A contour map of the maximum effective stress 
normalized with respect to the maximum pressure is shown in 
Fig. 5.2. The maximum effective stress is less than zero in 
the region directly beneath the load. Any cracks located in 
this area will lock. All cracks on the surface will also 
lock. The maximum effective stress increases away from the 
load and away from the surface. The highest fracture 
probabilities can be expected to occur near the locations of 
largest maximum effective stress. 

The probability of failure as a function of the maximum 
contact pressure is shown in Fig. 5.3. Bounds on the 
probability of failure are shown for three different 
correlation coefficients, p = 0, p = 1 and p = p(r..). Vhen 

J 
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0.002 

0.000 


Fig. 5.2 Contours of the maximum effective stress for 
normal loading only. 



0 



P 0 / E 


Fig. 5.3 Bounds on the probability of failure as a 

function of the normalized maximum pressure for 
different correlation coefficients. 



93 


the elements are not correlated, p = 0, the weakest link 
probability lies within the Ditlevsen bounds. However, when 
the system is fully correlated, the probability of failure 
is substantly lower than the weakest link probability, for a 
given load. This would indicate that for the whole system 
to fail more than one element needs to fail. An 
intermediate correlation function was assumed: 





5.4 


where r.. is the distance between the centroids of two 

J 

elements i and j and c is a constant, usually the average 
mesh size. This correlation is the most realistic because p 
is equal to one for two adjacent elements and decreases as 
the distance between elements increases. For this example c 
was assumed to be equal to one quarter of the width of the 
contact area or a/2. 

The purpose of establishing bounds on the probability 
of failure is to bracket a narrow range of expected failure 
probabilities for a given load. Reliability bounds are 
consistently used for design purposes where « 1. As the 
difference between bounds increases, generally with 
increasing load, they are not useful for a realistic 
bounding of the probability of failure. For example, for 
p - 0, as the probability of failure exceeds 0.5, the 
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difference between the bounds increase sharply and the 
concept of probability of failure is no longer important 
because failure has occurred. 

The largest failure probabilities occur in the area in 
which the maximum effective stress is greatest as shown by 
Fig. 5.4, a contour map of failure probabilities for each 
element. The normalized maximum pressure was 0.0165 when 
the map was drawn. A crack is most likely to propagate 
initially in the region of highest maximum effective stress. 
In practice, cracks originate beneath the surface and then 
grow around the inside of the bearing parallel to the 
surface and ultimately, the material peels off. This 
phenomenon is known as the shell effect. 

In the case of roller bearings, a frictional force is 
acting on the surface in the direction opposite to the 
relative motion. The tangential load is represented as a 
fraction of the pressure or 


q = fp 5.5 

where f is the friction coefficient and p is the pressure. 

A contour map of the maximum effective stress normalized 
with respect to the maximum pressure for a friction 
coefficient equal to 0.3 is shown in Fig. 5.5. 

As was the case when no friction is present cracks will 
lock in the area directly below the load. On the surface 
near the trailing edge of the tangential load the normal 



0.50 

0.10 


Fig. 5.4 Contour of the failure probabilities for each 
element . 
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0.010 

0.001 


Fig. 5.5 Contours of the maximum effective stress for 
normal and tangential loading. 
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stresses are tensile. The maximum effective stress in this 
area is the largest for the entire map. This dominate area 
will serve as the most likely place for the propagation of 
existing cracks. Another maximum on the map also occurs at 
the leading edge of the load. However ; it is not large 
enough to suspect that failure will occur at that point. 

The probability of failure as a function of the 
normalized maximum pressure is shown in Fig. 5.6. The 
probability of failure for the entire cylinder is 
approximately equal to the probability of failure of the 
element where the maximum effective stress is the highest. 
The upper and lower bound will converge to that value and 
are not dependent on the correlation coefficient. Fig. 5.7 
shows a contour map of the failure probabilities for each 
element. A comparison of Figs. 5.5 and 5.7 shows that the 
region of the highest maximum effective stress is also the 
region where the probability of failure is the highest. 

The probability of failure as a function of maximum 
normalized pressure is shown in Fig. 5.8 for different 
friction coefficients. As the friction coefficient 
increases the tensile stresses near the trailing edge become 
dominant and the failure probability increases dramatically 
for a given load. 



0.0030 


0.0035 


0.0040 


Fig. 5 


p 0 / e 


6 The probability of failure as a function of the 
maximum normalized pressure with normal and 
tangential loading. 
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Fig. 5.7 Contour of the failure probabilities for each 
element . 
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B. Notched Beam 

A test conducted at NASA Lewis Research Center whose 
purpose was to initiate a pre -crack for subsequent studies 
on fracture toughness involved applying a compressive load 
to a beam with a machined notch. The material chosen was a 
967. alumina ceramic whose failure data was available in the 
literature and is given in Table 5.1. The failure 
stresses found as a result of this analysis are not those 
expected for this experiment because the sample was fatigued 
for approximately 10 5 cycles before failure occurred. The 
purposes of this study were: to determine the location of 
failure, to define a loading limit under which failure will 
occur and to suggest a test configuration wherein a single 
crack grows under compressive loading at the base of the 
notch. 

A schematic showing the beams shape and loading 
condition is given in Fig. 5.9, where its height, 
h = 2.5 cm, width, w = 5.0 cm. and the thickness, 
t = 1.0 cm. The test configuration changed during the 
course of the experiment. The variable quantities were: the 
length of the notch, a, the length over which the load is 
applied, 1, and the distance between the back edge of the 
sample and the point of applied load, d. 

Initially, the material was loaded across the entire 
surface or d = 0, 1 = w and a = 2.7 cm. A crack was 
expected to grow at the base of the notch. Fig. 5.10 is a 
photograph of the beam at failure. After initial small 
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Fig. 5.9 Schematic of the loads applied to the notched 
beam . 
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Fig. 5.10 Tensile crack in the notched beam loaded across 
its entire width. 
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crack growth at the base of the notch, a second dominant 
crack originates at the top of the beam. The beam fails in 
tension. Vhen loading the specimen over three -fourths of 
its width, ultimate failure still occurs in the form of a 
tensile crack as shown in Fig. 5.11. 

A failure analysis on the beam was conducted by the 
author in order to determine if the mode of failure was 
justifyable. The finite element mesh is shown in Fig. 5.12. 
Because the beam is symmetric about its central axis, only 
half of the sample was modeled, however the entire beam is 
considered in the reliability analysis. The finite element 
model consists of 63 eight -node quadrilateral elements 
(plane -strain) , 230 nodes and 442 nodal degrees of freedom. 

After the displacements and stresses were analyzed 
using ADINA, the stresses and volumes are extracted from the 
output file and the failure for each element is determined. 
Vhen the weakest link probability of failure is equal to 
0.99, a contour map of the element probabilities of failure 
is drawn. The map for the beam loaded over its entire width 
is shown in Fig. 5.13. The element with the highest failure 
probability is located at the notch tip. However, almost 
directly above this element at the top of the beam, the 
region of highest tensile stress is located. The 
probability of failure in this region is significant in 
comparison with the failure probability at the notch tip. 

This explains the results of the actual experiment. 

The crack growth initiating at the notch tip is shown by the 
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Fig. 5.11 Tensile crack in the notched beam loaded across 
three quarters of its entire width. 





x notch tip 


Fig. 5.12 Discretization of the notched beam. 
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Fig. 5.13 Contour map of the probabilities of failure for 
the beam loaded over its entire width. 
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highest probability of failure occurring at that location. 
Due to the nature of brittle fracture in compression, the 
crack growth was local. As additional load is applied, the 
regions of tensile stress will dominate. The tensile crack 
which breaks the sample is initiated at the top. 

Another test was conducted under the following 
conditions: a = 2.7 cm, 1 = 1.12 cm, and w = 1.4 cm. The 
results of the experiment show that a dominant crack 
initiates in the region of tensile stress. A photograph of 
that crack is shown in Fig. 5.14. A small crack is 
initiated at the base of the notch in compression with 
ultimate failure caused by the tensile crack. 

Vhen the weakest link probability of failure is equal 
to 0.99, a contour map of the failure probabilities for 
individual elements is shown in Fig. 5.15. The probability 
of failure is the largest for the element at the notch tip, 
however at the top of the beam away from the notch, the 
probability of failure is significant for the tensile stress 
region. 

Fig. 5.16 shows the probability of failure as a 
function of the applied load. The weakest link probability 
is shown along with the bounds for three correlation 
coefficients, p = 0, p = 1 and p = p(r^) as given by 
eq. 5.4 with c = 0.3 cm. Vhen the elements are not 
correlated, p = 0, the weakest link probability is 
approximately equal to the Ditlevsen bound, yielding a 
unique solution. However, when the system is fully 
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correlated, the difference between the bounds is substantly 
greater than for an uncorrelated system. The spread between 
the bounds for the intermediate correlation coefficient lies 
between that for a fully correlated and uncorrelated system 
for a given load. These trends are the same as those found 
for the cylinders in contact as given by Fig. 5.3. 

Vith an applied load of 44.5 kN (10,000 lbs) the 
tensile crack appears after 97,320 cycles. This load is 
approximately 807. of the load required to yield a 
probability of failure of 0.01 for the weakest link model as 
shown in Fig. 5.16. 

In order to determine a test configuration wherein a 
crack grows at the base of the notch, a sensitivity study of 
the possible loading configurations was completed. Vith the 
beams height to width ratio of 0.5, an improved test 
condition was not found. However, if the height was 
doubled, the tensile stresses away from the notch decrease. 
Then the probability of failure in these regions will 
decrease. A contour map of the probabilities of failure of 
a beam loaded as follows: a = 2.7 cm, 1 = 1.12 cm, and 
w = 1.4 cm is shown in Fig. 5.17. The height of this 
specimen is doubled, h = 5.0 cm. The probability of failure 
in the elements away from the notch tip has decreased. 

The probability of failure as a function of the applied 
load for the two different heights (identical loading 
conditions), is shown in Fig. 5.18. For a given load the 
failure probability decreases for the beam with h = 5.0 cm. 


PRECEDING PAGE blank not filmed 
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Fig. 5.17 Contour map of the probabilities of failure for 
the beam loaded over a section of its entire 
width where the height is doubled. 
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Fig. 5.18 The weakest link probability of failure as a 
function of the applied load for different 
beam heights. 
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The tensile stress in the beam decreases. The volume or the 
surface area increases for the beam vhose height is greater. 
The combined effect of the two results is the lowering of 
the probability of failure at a given load for the beam 
which is larger. 

In conclusion, after analyzing the probability of 
failure for different loading conditions and beam heights, 
the configuration which would most likely result in crack 
growth at the notch tip and not in tension is given in 
Table 5.4. Since the probability of failure decreased for 
this specimen at a given load, the load which was applied in 
the earlier experiments may be increased slightly or the 
number of cycles needed for compressive failure at the notch 
tip may be increased without the risk of remote tensile 
failure. 


Table 5.4 Dimensions and loading condition of the 
notched beam. 


Dimensions 


width = 5.0 cm (1.92 in.) 

height = 5.0 cm (1.92 in.) 

thickness = 1.0 cm (0.40 in.) 

a = 2.7 cm (1.064 in.) 

Loading condition 

d = 1.12 cm (0.44 in.) 

1 = 1.40 cm (0.55 in.) 

Applied load = 44.5 kN (10,000 lbs) 



CHAPTER VI 


CONCLUSIONS 


A . Summary 

The objective of this study vas to analyze the 
reliability of brittle materials under contact stress 
conditions. The analysis has shown the following 
conclusions : 

- The failure of brittle materials whose compressive 
strength is much larger than its tensile strength, is 
analyzed using a Batdorf model modified to include the 
reduction in shear due to the effect of the compressive 
stresses on the crack face. 

- Fracture of brittle materials is modeled as a series 
system. 

- As the system becomes more fully correlated the 
probability of failure decreases for a given load. 

- The spread of the bounds increases as the correlation 
coefficient increases. 

- The nature of brittle crack growth in compression 
results in local crack growth. The presence of a similar 
crack in tension will lead to catastrophic failure. Failure 
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of brittle materials is biased in tension. 

- The probability of failure is largest in the regions 
where the maximum effective stress is greatest. 

B . Further Vork 

- A reliability analysis including both the effect of 
the shear and crushing compressive crack growth mechanisms 
should be considered. 

- An experimental study should be conducted to 
determine the internal friction coefficient as it relates to 
characterizing the crack density function for brittle 
materials . 

- Nonlinear fracture envelopes should be applied in 
examining the reliability of brittle materials under 
compressive loading. 
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